Solving the Schroedinger equation for the Sherrington-Kirkpatrick model in a 

transverse field 
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By numerically solving the Schroedinger equation for small sizes we investigate the quantum critical 
point of the infinite-range Ising spin glass in a transverse field at zero temperature. Despite its 
simplicity the method yields accurate information on the value of the critical field and critical 
exponents. We obtain r c = 1.47 ± 0.01 and check that exponents are in agreement with analytical 
approaches. 
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There has recently been renewed interest in the study 
of quantum phase transitions in disordered systems [jy . 
In particular, Ising spin glass models in a transverse field 
are simple systems in which to study the effect of com- 
petition between randomness and quantum fluctuations. 
The case of infinite-range models is especially interesting 
because they show non-trivial quantum phase transitions 
yet are to some extent amenable to analytical computa- 
tions. The canonical example in this family of models 
is the quantum Sherrington-Kirkpatrick (SK) model in a 
transverse field. At zero transverse field this reduces to 
the usual SK model which has a finite temperature transi- 
tion to low temperature phase where replica symmetry is 
broken || . As the transverse field is turned on, spin-glass 
ordering occurs at lower temperatures and above a cer- 
tain critical field the spin-glass order is completely sup- 
pressed at the expense of ordering in the transverse direc- 
tion. Our understanding of this model was significantly 
extended by the non-perturbative analysis of Miller and 
Huse H and also by the different approach of Ye, Sachdev 
and Read The critical behavior is now well estab- 
lished, and values for exponents, predictions for logarith- 
mic corrections and estimates of the value of the critical 
field are known. The model is therefore well adapted 
as a testing ground for numerical methods to investigate 
quantum phase transitions. From this point of view, the 
phase diagram of the quantum SK model in a transverse 
field, and its zero temperature critical behavior have been 
studied using numerical techniques such as spin summa- 
tion [^) , perturbation expansions || and quantum Monte 
Carlo methods ||. It is the purpose of this letter to in- 
troduce a new numerical approach based on the intuitive 
method of directly solving the Schroedinger equation for 
finite systems. Despite its simplicity, this method is able 
to give quantitative information on the value of the crit- 
ical field and critical exponents even for the very small 



size systems we consider. 

The SK model in a transverse field is defined by the 
Hamiltonian, 



H = H - TM X 



(i) 



where of , af are Pauli spin matrices and T is the the 
transverse field. The indices i, j run from 1 to N where N 
is the number of sites. The Jy are Gaussian distributed 
random variables with zero mean and 1/N variance. Ho 
is the term we call the interaction energy, while M. x 
stands for the magnetisation in the transverse direction. 

We propose to study this model by the direct method 
of numerically solving the real time Schroedinger dynam- 
ics, 



-win 



(2) 



The wave function, \tp(t)}, of the system at time t can be 
written as a linear combination of basis states, 



w*)> = 5>.(*)i"> 



(3) 



We have chosen the basis states, {\v); v = 1, .., 2 N } to be 
eigenstates of each of the spin operators — 1,-^V}, 

\v) = \si, S2 ■ ■ ■ snY- This choice gives a geometric mean- 
ing to equation (||) because these eigenstates can also 
be interpreted as the vertices of a unit hypercubic cell 
of dimension TV. Each vertex of this hypercubic Hilbcrt 
space is assigned a label v and a corresponding complex 
variable a V: which together define the state of the sys- 
tem. This geometric picture can also be used to under- 
stand the action of the Hamiltonian operator on the basis 
states \v). The action of the first term in equation (111) on 
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the state \v) is diagonal with eigenvalue E® which is pre- 
cisely the energy of the classical SK model in that state 
(E® — {is\Ha\v))- The operator a l x acting on a given 
eigenstate \v) changes the value of one spin which cor- 
responds to an adjacent vertex of the hypercube. The 
dynamical equations for the a v become, 

i^ = E°a v (t)-TY,a»(t) (4) 
/*(«/) 

where n(v) are nearest neighbours to the vertex v in the 
hypercubic cell. The geometric picture also facilitates 
efficient computer code for this problem. 

We wish to calculate thermodynamic properties of the 
Hamiltonian ([!]) at zero temperature. Such information 
could be obtained by finding the static ground state of the 
Hamiltonian H. However, because we want this ground 
state for a range of transverse fields, it is convenient to 
use a dynamical procedure. At large T the Hamiltonian 
can be diagonalised and the ground state is given by: 
a u (t = 0) = 1/2 3". Starting from this configuration we 
reduce the transverse field adiabatically slowly, thus en- 
suring that the system remains in its ground state. This 
procedure is recommended by its simple physical inter- 
pretation but is not the most efficient method that could 
be envisioned. For example the phase of the wavefunc- 
tion is not of interest, and for large systems a gain in 
speed could be achieved by some less direct method. 

We have numerically integrated eq. ([|) for different val- 
ues of N ranging from N = 2 to N = 13, using a simple 
Euler algorithm. The value of the time step can be fixed 
by testing the conservation of energy for some excited 
state at fixed T. For the ground state we can be less 
careful and we choose a time step dt = 0.01. The trans- 
verse field is allowed to decrease linearly from T = 3 
down to r = 0. We find that a total time of 100 units 
(amounting to 10000 integration steps) gives sufficiently 
slow variation of T for the adiabatic theorem to hold. 
The method has also been checked against the analytic 
solution of the model for N = 2. The errors from the 
discretisation of the Schroedinger equation, from the adi- 
abatic approximation and from the finite initial value of 
r are therefore well under control. The main source of 
error comes from sample to sample fluctuations. Data 
was averaged over many samples; ranging from 50000 for 
the smallest systems (N = 3, 4, 5) to 3000 for the largest 
sizes (N — 10, 11, 13). We have also considered 100 sam- 
ples at N = 17 to confirm the tendency of the data, but 
have not used these points in our fits due to the large 
errors. 
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FIG. 1. Interaction energy against tranverse field for (from 
top to bottom) N = 2 (analytic curve), N — 5 (50000 sam- 
ples), N = 8 (30000 samples) and N = 13 (3000 samples). 
Errors are not shown, but are always less than 10~ 3 . The low- 
est continuous line is the extrapolated data, Eo(oo), and the 
points have been obtained by quantum Monte Carlo methods 
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FIG. 2. Interaction Energy against 1/N. At F — 2.5 in 
the QP phase (top) and F = 1.45 near the QC point (bottom). 
The points have error bars and the continuous curves show 
the power law fits. 



The simplest variable is the interaction energy, Eq = 
(0\Ho\0) = J2 U E v a v a*. It is plotted in figure 1 as a func- 
tion of r for several different sizes. In figure 2 we show 
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Eq for three different values of T as a function of 1/N. 
The different values of T are T = 2.5 in the quantum 
paramagnetic phase (QP) and T = 1.45 near the quan- 
tum critical (QC) point (see below). Data has been fitted 
with the least squares method to a power law of the type 
E (N) = E Q (oo)+aN- b |,|. In the QP phase b ~ 0.93, 
so corrections are essentially 1/N as expected. Close to 
the QC point we find b = 0.73 ± 0.02 in agreement with 
the expected mean- field exponent b = 1 + -f- + = 3/4 
where v is the correlation length exponent (v = 1/4), z 
is the dynamical exponent (z = 2) and d u is the upper 
critical dimension (d u =8). It is remarkable that even 
very small size systems fit on the curve. This data is 
summarized in figure 1 where we have also shown the 
best fit parameters Eq{oo) as a function of T in the re- 
gion r > 1.2. The points are numerical data obtained by 
independent Quantum Monte Carlo calculations |)| and 
show reasonable agreement with the extrapolated values. 
At smaller T, b decreases and it is no longer possible to 
ignore sub-leading corrections; nonetheless, at T = we 
find E {oo) ~ —0.763 in agreement with the theory 

These results give us confidence in the method and 
encourage us to investigate the transition more closely. 
Because of the spin glass nature of the transition, a clear 
signal does not appear in the more ordinary thermody- 
namic functions. A divergence will occur in the non- 
linear susceptibility and we have also seen a peak in a 
Binder-like parameter for the kurtosis of the sample to 
sample distribution of the interaction energy. These ob- 
servations suggest a transition in the region of T ~ 1.5 
as expected, but are not the best way to obtain accurate 
information. To determine the critical field and criti- 
cal exponents we consider the longitudinal susceptibility 
associated to the magnetisation Ai z = of. It can be 
shown jl(]] that the longitudinal susceptibility (x) for the 
SK model is precisely equal to 1 at the QP-QG boundary. 
To this end we have numerically computed \. The usual 
quantum mechanical formula, 

* = E« (5) 

where \n) denotes the energy eigenstate with energy E ni 
is not adapted to our needs because it requires knowl- 
edge of the full spectrum. Instead, we solve the problem 
directly by applying a longitudinal magnetic field h small 
enough to be in the linear response regime. The suscepti- 
bility is computed as x — {0'\M z \0')/h where |0') stands 
for the ground state in the presence of the field h. We 
need only solve the Schroedinger equation once for the 
perturbed Hamiltonian TL + hAi z , since the magnetisa- 
tion of the unperturbed problem (without magnetic field) 
is strictly zero, as follows from a quantum mechanical 
symmetry. For zero transverse field this symmetry is the 
spin reversal symmetry of the classical SK model. In the 



thermodynamic limit N — > oo only one state is selected, 
but for our finite systems the wavefunction always con- 
tains a mixture of opposite magnetisation states. In the 
perturbed case this symmetry is lost and at small trans- 
verse fields the action of the perturbation is to shift the 
wavefunction to a single magnetised state. This quantum 
tunneling introduces a new time scale into the problem 
and the parameters we use in solving the Schroedinger 
equation should be re-examined. We have checked that 
the parameters we use give small errors (less than the er- 
rors from sample fluctuations) around the region of crit- 
icality. 
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FIG. 3. Critical Transverse field against 1/N. 

The actual value of the perturbing field h can be taken 
in a wide range without affecting results, and we present 
data for h as small as 10~ 7 . Using the exact condition 
for criticality we define the critical transverse field by 
%(r) = 1, and in figure 3 we plot this T against 1/N. 
Fitting the data to a power law behavior of the type 
r = T c + aN- b we find T c = 1.47 ± 0.02, a = -0.485 ± 
0.002, b = 0.53 ± 0.02 in good agreement with the results 
obtained in perturbation theory by Ishii and Yamamoto 
Q (r c = 1.506) and with the result obtained by Miller 
and Huse || (r c = 1.46 ± 0.01). The coefficient b is very 
close to the expected value b = -4- = 1/2. If we consider 
the scaling at T = 1.47 we find that \ converges to 1 as 
N~ c with a value of the exponent c ~ 0.29 compatible 
with c = j- = 1/4. The data collapse in the scaling 

U 2 

region is shown in figure 4 where (1 — \)N du is plotted 
as a function of N~^(T — T c ). The collapse is good 
and confirms the expected values v = 1/4 and d u = 8. 
Even the result for N — 2 lies close to the collapse line. 
Considering the simplicity of the method and the small 
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sizes considered the matching with data reported in the 
literature is impressive. This is particularly the case since 
logarithmic corrections are known to be present |]J4|]. 



Z 0.3 - 



[10] 

[11] 
[12] 



1987); K.H. Fischer and J. A. Hertz, Spin Glasses (Cam- 
bridge University Press, Cambridge 1991). 
J. Miller and D. Huse, Phys. Rev. Lett. 70, 3147 (1993). 
J. Ye, S. Sachdev and N. Read, Phys. Rev. Lett. 70, 4011 
(1993). 

G. Biittner and K.D. Usadel, Phys. Rev. B 41, 428 (1990); 
Y.Y. Goldschmidt and P.Y. Lai, Phys. Rev. Lett. 64, 
2467 (1990). 

H. Ishii and Y. Yamamoto, J. Phys. C 18, 6225 (1985); 
J. Phys. C 20, 6053 (1987). 

G. Parisi, F. Ritort and F. Slanina, J. Phys. A (Math. 

Gen.) 26 (1993) 247; 26 (1993) 3775. 

S. Kirkpatrick and A.P. Young, Phys. Rev. B 25 (1982) 

440. 

J.V. Alv arez and F. Ritort, Preprint |coiid-mat 
9603012] to appear in J. Phys. A. 

A.J. Bray and M.A . Moore, J. Phys C. 1 3, L655 (1980). 



F. Ritort, Preprint |cond-mat 9607044 . 

T. Senhil and S.N. Majumdar, Phys. Rev. Lett. 76, 
(1996). 



3001 



N 1/a (r- r c ) ' 

FIG. 4. Data collapse for N = 5 (triangles), N = 8 
(squares) and N = 13 (circles). Shown for the range 
1.2 < T < 2.0 using our value of F c = 1.47. The continu- 
ous line is the analytical result for N = 2. 



In summary, a new and simple numerical method yield- 
ing good estimates of the critical field and confirming the 
critical exponents for the quantum phase transition of 
the SK model has been reported. The method consists 
in solving the Schroedinger equation for small sizes and 
computing expectation values in the ground state of the 
system. Quite remarkably, the system is within the scal- 
ing region even for very small sizes. The application of 
this method to other disordered systems such as the ran- 
dom orthogonal model |ll| and the quantum Potts glass 
Jl2| would be very welcome. 
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